# Daten werden je nach Sprache zentriert (silent reading speed)
data <- data %>%
  mutate(words_per_minute_centered = case_when(
    language == "German"  ~ words_per_minute - 260,
    language == "English" ~ words_per_minute - 238
  ))



############################################################################## #
## With random intercept                                                    ####
############################################################################## #

fit_gamma_RE <- brm(
  words_per_minute ~ 1 + (1 | study),
  data = data,
  family = Gamma(link = "log"),
  chains = 4, iter = 4000, cores = 4, seed = 123
)

# --- get posterior mean and CI ---
# Posterior: erwarteter Populationsmittelwert (auf Originalskala)
newdata <- data.frame(study = NA)   # NA in Gruppierungsvariable => neue Ebene (re_formula = NA fixiert random effects auf 0)
epred <- posterior_epred(fit_gamma_RE, 
                         newdata = newdata, # predict these cases
                         re_formula = NA)  # ignore random effects, gets population average (only fixed effects)
epred_vec <- as.numeric(epred[,1])         # get vector with posterior draws

mean(epred_vec)                       # Posterior mean of population expected reading_speed
quantile(epred_vec, c(0.025, 0.975))  # 95% CI

# --- hypothesis() on log-Scale ---
log_500 <- log(500 - 238) # words per minute for skim reading and reading limit 
log_900 <- log(900 - 238) # aren't clearly defined. We therefore use the English
                          # language value (and not both languages separately)
                          # to subtract from 500 and 900

hypothesis(fit_gamma_RE, "Intercept > 0")
hypothesis(fit_gamma_RE, paste0("Intercept > ", log_500))
hypothesis(fit_gamma_RE, paste0("Intercept > ", log_900))





############################################################################## #
## Without random intercept                                                 ####
############################################################################## #
library(brms)

# Gamma-Modell mit log-Link, nur Intercept
fit_gamma_noRE <- brm(
  words_per_minute ~ 1,
  data = tmp,
  family = Gamma(link = "log"),
  chains = 4, iter = 4000, cores = 4, seed = 123
)


# --- get posterior mean and CI ---
epred <- posterior_epred(fit_gamma_noRE)
epred_vec <- as.numeric(epred[, 1])

mean(epred_vec)                       # Posterior mean
quantile(epred_vec, c(0.025, 0.975))  # 95%-Intervall

# --- hypothesis() on log-Scale ---
log_500 <- log(500 - 238) # words per minute for skim reading and reading limit 
log_900 <- log(900 - 238) # aren't clearly defined. We therefore use the English
                          # language value (and not both languages separately)
                          # to subtract from 500 and 900

hypothesis(fit_gamma_RE, "Intercept > 0")
hypothesis(fit_gamma_RE, paste0("Intercept > ", log_500))
hypothesis(fit_gamma_RE, paste0("Intercept > ", log_900))

library(brms)

fit_mediation <- brm(
  data = mydata,
  formula = 
    bf(dwell_time ~ experience_participant + experience_conducting + trust_benevolence + perceived_risk,
       family = Gamma(link = "log")) +  # log link as data are surely skewed
    bf(reading_engagement ~ experience_participant + experience_conducting + trust_benevolence + perceived_risk,
       family = student()) +
    bf(understanding ~ experience_participant + experience_conducting + trust_benevolence + perceived_risk +
         dwell_time + reading_engagement,
       family = student()) +
    set_rescor(TRUE),
  chains = 4, iter = 4000, cores = 4, control = list(adapt_delta = 0.95)
)

posterior <- as_draws_df(fit_mediation)

# Indirect effects over dwell time ############################################
indirect_dwell <- posterior$b_dwelltime_experience_participant * posterior$b_understanding_dwell_time +
                  posterior$b_dwelltime_experience_conducting * posterior$b_understanding_dwell_time +
                  posterior$b_dwelltime_trust_benevolence * posterior$b_understanding_dwell_time +
                  posterior$b_dwelltime_perceived_risk * posterior$b_understanding_dwell_time

# Credible Interval + Posterior mean
quantile(indirect_dwell, c(0.025, 0.975))
mean(indirect_dwell)


# Indirect effects over reading engagement #####################################
indirect_engagement <- posterior$b_readingengagement_experience_participant * posterior$b_understanding_reading_engagement +
                  posterior$b_readingengagement_experience_conducting * posterior$b_understanding_reading_engagement +
                  posterior$b_readingengagement_trust_benevolence * posterior$b_understanding_reading_engagement +
                  posterior$b_readingengagement_perceived_risk * posterior$b_understanding_reading_engagement

# Credible Interval + Posterior mean
quantile(indirect_engagement, c(0.025, 0.975))
mean(indirect_engagement)


# total indirect effects #######################################################
indirect_total <- indirect_dwell + indirect_engagement

# Credible Interval + Posterior mean
quantile(indirect_total, c(0.025, 0.975))
mean(indirect_total)

library(brms)

fit <- brm(
  formula = bf(willingness_participate ~ data_sharing * data_type_interview + 
                                         data_sharing * data_type_video +
                                         data_sharing * perceived_sensitivity),
  data = mydata,
  family = gaussian(),
  chains = 4,
  cores = 4,
  iter = 4000,
  warmup = 1000,
  control = list(adapt_delta = 0.95))